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Abstract. We present constraints on neutrino masses, the primordial fluctuation spectrum from 
inflation, and other parameters of the ACDM model, using the one-dimensional Fycr-forest power 
spectrum measured by Palanque-Delabrouille et al. [1] from the Baryon Oscillation Spectroscopic 
Survey (BOSS) of the Sloan Digital Sky Survey (SDSS-III), complemented by Planck 2015 cosmic 
microwave background (CMB) data and other cosmological probes. This paper improves on the 
previous analysis by Palanque-Delabrouille et al. [2] by using a more powerful set of calibrating 
hydrodynamical simulations that reduces uncertainties associated with resolution and box size, by 
adopting a more flexible set of nuisance parameters for describing the evolution of the intergalactic 
medium, by including additional freedom to account for systematic uncertainties, and by using Planck 
2015 constraints in place of Planck 2013. 

Fitting Fycr data alone leads to cosmological parameters in excellent agreement with the values 
derived independently from CMB data, except for a weak tension on the scalar index Combining 
BOSS Fycr with Planck CMB constrains the sum of neutrino masses to ^ niy < 0.12 eV (95% C.F.) 
including all identified systematic uncertainties, tighter than our previous limit (0.15 eV) and more 
robust. Adding Fycr data to CMB data reduces the uncertainties on the optical depth to reionization 
T, through the correlation of r with cr%. Similarly, correlations between cosmological parameters 
help in constraining the tensor-to-scalar ratio of primordial fluctuations r. The tension on can be 
accommodated by allowing for a running dn^/dlnk. Allowing running as a free parameter in the fits 
does not change the limit on ^ niy. We discuss possible interpretations of these results in the context 
of slow-roll inflation. 
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1 Introduction 

The flux power spectrum of the Lyman-a (Lya) forest in quasar absorption spectra is a powerful 
tool to study clustering in the Universe, at redshifts ~ 2 - 4. Compared to a model derived from 
a set of dedicated hydrodynamical simulations, the Lya-flux power spectrum can provide valuable 
information on the formation of structures and their evolution. Furthermore, by probing scales down 
to a few Mpc, the ID flux power spectrum is also sensitive to neutrino masses through the suppression 
of power on small scales that neutrinos induce because they become non-relativistic at small redshift 
and they therefore free-stream during most of the history of structure formation. We here use the 
ID Lya flux power spectrum measured with the DR9 release of BOSS quasar data, and a grid of 36 
hydrodynamical simulations having a resolution equivalent to 3 x 3072^ particles in a (100 Mpc)^ 
box, to constrain both cosmology and the sum of the neutrino masses ^ ruy. 

Cosmic Microwave Background (CMB) can also constrain ^ my. In the standard thermal his¬ 
tory of the Universe, massless neutrinos have an average momentum corresponding to 3.15Ty ~ 
0.58 eV at the epoch of last scattering. For X > 3 x 0.58 = 1.7 eV, the neutrinos are still relativis¬ 
tic at recombination, and have no significant impact on the primary CMB anisotropies. However, for 
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any mass value, neutrinos leave a signature on the CMB angular power spectrum through the inte¬ 
grated Sachs-Wolf effect and through lensing [3, 4]. The latest limit on ^ nty from CMB data alone 
is at the level of 0.7 eV [5]. 

Lya data alone have sensitivity to X at the level of about 1 eV due to the fact that the 
scales probed by Lya forests are in the region where the ratio of the power spectra for massive to 
massless neutrinos is quite flat (cf. Figure 1). However, a tight constraint on ^ ruy can be obtained 
by combining CMB data, which probe the initial power spectrum unaffected by ^ ruy, and Lya data, 
which probe the suppressed power spectrum. Thus, Lya measures the power spectrum level, defined 
by cs and Q-m, CMB provides the correlations between these parameters and Yj ^v, and the joint 
use of these two probes significantly improves the constraint on Y compared to what either probe 
alone can achieve. 



Figure 1: Linear theory prediction for the matter power spectra with massive neutrinos, normalized 
to the corresponding massless neutrino case. The grey zone delimits the range ofk covered by the ID 
Lya flux power spectrum from the BOSS survey. 


The layout of the paper is as follows. The first part of section 2 presents the upgrades in the Lya, 
CMB and Baryon Acoustic Oscillation (BAO) data sets used for this work. The second part summa¬ 
rizes a number of improvements in the methodology: changes in the accounting of the uncertainties 
of the hydrodynamical simulations, and updates of the likelihood parameters to allow for additional 
freedom in the IGM model or in the instrumental systematic effects. The main objective of section 3.1 
is to present what Lya data alone have to say about cosmology. The base model we consider is a flat 
ACDM cosmology with massive neutrinos, thereafter referred to as the base ACDMv cosmology. 
We start by giving the constraints measured on the five relevant parameters (eg, Us, Llm, Hq, Y ^v), 
and we briefly discuss the values of the ‘nuisance’ parameters. In section 3.2, we include additional 
data, namely several configurations of CMB data and, occasionally, BAO measurements. We present 
the results obtained on the parameters of our base ACDMv cosmology with various combinations of 
these data sets. Finally, we discuss extensions to the base ACDMv cosmology. We present how Lya 
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data can contribute to constraining additional parameters through their correlations with parameters 
that Lya data are sensitive to; we thus give constraints on the reionization optical depth in section 3.3, 
and on primordial fluctuations (e.g., the running of the scalar spectral index dn^/dlnA: and the ratio 
of tensor to scalar modes r) in section 3.4. We discuss the small impact of the running of on the 
constraints on Yj We do not discuss here the combined Lyor + CMB constraint on the number of 
neutrino species N^ff, since this is the object of a dedicated study [6]. 

This paper refers extensively to the earlier paper [2] that reported the first constraints on cos¬ 
mological parameters and total neutrino mass using Lycr data from the SDSS-III/BOSS survey. To 
simplify the presentation and make it easier for the reader to identify any reference to this earlier pa¬ 
per, we will henceforth refer to it as Paper I. We also refer the reader to [7] for a detailed description 
of the grid of hydrodynamical simulations used in this work, and to [8] for the implementation of 
neutrinos and their impact on the ID flux power spectrum. Definitions of the most relevant symbols 
used in this paper can be found in Tables 1 and 2. 

Table 1: Definition of astrophysical parameters 


Parameter Definition 

6 = pi (p) . Normalized baryonic density p of IGM 

T . Temperature of IGM modeled by T = Tq- SP'~^ 

Tq . Normalization temperature of IGM at z = 3 

y . Logarithmic slope of 6 dependence of IGM temperature at z = 3 

77^0 . Logarithmic slope of redshift dependence of Tq (different for z < or > 3) 

rp' . Logarithmic slope of redshift dependence of y 

A’’. Effective optical depth of Lya absorption at z - 3 

if . Logarithmic slope of redshift dependence of A^ 

/siiii. Fraction of Si III absorption relative to Lya absorption 

/siii . Fraction of Si II absorption relative to Fya absorption 


2 Data and methodology 

The methodology, notation and parameters used in this paper mostly follow those described in Pa¬ 
per I [2]. We, however, include additional data, as described in Sec. 2.1. We also present a number 
of updates in the methodology. We include additional parameters in the likelihood to account for all 
the systematic effects we identified, whefher relafed fo fhe insfrumenf or fo the simulations. We have 
a revised estimate of the sampling variance, and an improved model of the impact of the splicing 
technique on the simulated power spectrum. We also allow for additional freedom in the model of 
the intergalactic medium (IGM) temperature. All these changes are described in Sec. 2.2. 

2.1 Changes in the data 

This paper is mostly focused on the analysis of Fya data to provide constraints on cosmological pa¬ 
rameters. We also include at times additional probes, namely measurements of the cosmic microwave 
background and of large-scale structures at redshifts z < 1, to improve the constraints we derive. 
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Table 2: Definition of cosmological parameters 


Parameter Definition 

Q.,„ . Matter fraction today (compared to critical density) 

Hq . Expansion rate today in km s“^ Mpc“' 

X nty. Sum of neutrino masses in eV 

CTg . RMS matter fluctuation amplitude today in linear theory 

T . Optical depth to reionization 

Zre . Redshift where reionization fraction is 50% 

As . Scalar power spectrum amplitude 

ris . Scalar spectral index 

ko = 0.05 Mpc“^ . Pivot scale of CMB 

dris/dlnk . Running of scalar spectral index 

r. Tensor-to-scalar power ratio 

£1 = -HjH^ . First Hubble hierarchy parameter 

£,'+1 =-e,7(//e,) . (i+l)st Hubble hierarchy parameter 


2.1.1 Lya 

As our large-scale structure probe, we use the ID Lya-flux power spectrum measurement [1] from 
the first release of BOSS quasar data. The data consist of a sample of 13 821 spectra selected from 
the larger sample of about 60 000 quasar spectra of the SDSS-III/BOSS DR9 [9-14] on the basis of 
their high quality, high signal-to-noise ratio and good spectral resolution (< 85 km s“^ on average 
over a quasar forest). In Paper I, we chose to focus the analysis on the first ten redshift bins, spanning 
the range 2.1 < z < 4.1. With the improvement on our model of the IGM temperature, however, 
we are now able to use all 12 redshift bins, thus extending the redshift coverage of the analysis to 

2.1 < z < 4.5 as shown on Fig. 2. We thus now do the analysis on 420 Fya data points, consisting of 
12 redshift bins and 35 k bins. 

2.1.2 Cosmic microwave background 

In Paper I, the cosmic microwave background (CMB) data and results we used were described in 
the March 2013 Planck cosmological parameters paper [15]. When combining with CMB data, we 
now use the more recent cosmological results obtained from the full Planck mission and presented 
in [5]. We consider several subsets of Planck data. The base configuration, denoted ‘TT-i-lowP’ 
as in [5], uses the TT spectra at low and high multipoles and the polarization information up to 
multipoles { - 29 (‘lowP’). We also use at times the configuration based on TT, TE and EE spectra, 
along with the low-multipole polarization, denoted ‘TT-i-TE-i-EE-i-lowP’. When studying inflation 
parameters, we will also use the BICEP2/keck Array-Planck joint analysis [16], which combines 
the high sensitivity B-mode maps from BICEP2 and Keck Array with the Planck maps at higher 
frequencies where dust emission dominates. This set will be denoted ‘BKP’. 

2.1.3 Baryon acoustic oscillations 

Finally, we occasionally combine CMB data with measurements of the BAO scale by 6dFGS [17], 
SDSS main galaxy sample [18], BOSS-EOWZ [19] and CMASS-DRll [19]. Theses measurements 
are henceforth globally denoted ’BAG’. The additional constraints that these measurement provide 
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Figure 2: ID Lya forest power spectrum from the SDSS-llI/BOSS DR9 data. The solid curves show 
the best-fit model obtained in section 3.1 when considering Lya data alone. The oscillations in 
these model predictions (and, presumably, in the measurements) arise principally from Lya-Si III 
correlations, which occur at a wavelength separation AA = 9.2A. 


on cosmological parameters are included in the present work with their full correlation with CMB 
data. Both CMB and BAO constraints are taken from the Markov Chains publicly available through 
the official Planck Legacy Archive at http://pla.esac.esa.int. 

2.2 Changes in the methodology and models 

We interpret the ID Lya-llux power spectrum using a likelihood built around three categories of 
parameters which are floated in the minimization procedure. The first category describes the cosmo¬ 
logical model in the simplest case of ACDM assuming a flat Universe. The second category models 
the astrophysics within the IGM, and the relationship between the gas temperature and its density. 
The purpose of the third category is to describe the imperfections of our measurement of the ID 
power spectrum. This likelihood allows us to compare the measurement to the power spectrum pre¬ 
dicted from hydrodynamical simulations. The changes in the simulation model or in the likelihood 
compared to Paper I are described below. 

2.2.1 Sample variance 

A sample variance is expected on large scales since the size of the simulation volume is similar to the 
largest modes measured. We improved our estimate of its contribution to the simulation uncertainties 
by computing, for each mode, the variance of the difference from average of the ID Lya-flux power 
spectrum for five simulations run with exactly the same cosmological and astrophysical parameters 
but different random seeds to initiate the distribution of particles (cf. Fig. 3). As expected, this 
test shows an excess of variance at small k, compared to the uncertainty measured within each run, 
which we model by a function of the form (a -v b exp{-ck))^ where k is the wavenumber, a - 0.004, 
b = 0.023 and c = -356.6. This additional variance is added in quadrature to the simulation statistical 
variance. 
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Figure 3: Visualization of the baryonic gas at z = 2.2 for five simulations run with identical param¬ 
eters but different random seeds to compute the initial conditions. Simulations are using 2 x 768^ 
particles in a (100 h^^ Mpcf box. Color represents gas temperature (from blue to red) and density is 
mapped to intensity. The plot at the bottom right illustrates the relative total uncertainty on the data 
power spectrum (grey) and the level of sampling variance (colored dots for each redshift bin, same 
color-code as Fig. 2). 
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2.2.2 Splicing 


The ID power spectrum is computed using a set of three simulations of different mass resolutions 
(192^ or 768^ particles per species) and sizes (boxes 25 or 100 h~^ Mpc on a side). A splicing 
technique [20] is then used to combine the power spectra and provide the equivalent of a power 
spectrum measured from a single simulation spanning a total volume of (100 h~^ Mpc)^ with 3072^ 
particles per species. We study the accuracy of the method by comparing a spliced power spectrum 
^spliced built from a simulation set with 128^ or 512^ particles per species in a box of side length 25 or 
100 Mpc, with the corresponding ‘exact’ power spectrum /"exact obtained from a simulation with 
2048^ particles per species in a box of side length 100 h~^ Mpc.^ The residuals /’spiiced/Z’exact are 
shown in figure 4. They are modeled by a broken line, with pivot point at the mode kp corresponding 
to the change of regime from constant to k-dependent resolution correction in the splicing technique 
(cf. details in [7]). The vertical offset at kp is let free, the slope at k < kp is fixed to -21 or -27 km.s“' 
for z < 3.5 or z > 3.5 respectively, and the slope at k > kp is let free. We allow for a small redshift 
dependence in the splicing correction (decreasing residual offset at the pivot point, and increasing 
slope in the k > kp region as z increases). This model is an improvement over Paper I where the 
accuracy of the splicing technique was estimated from less resolved simulations and modeled by a 
single redshift-independent linear function of k over all modes. 


N = 2048 



Figure 4: Residuals of the spliced to the exact power spectrum for z - 2.2, 3.4 and z = 4.4 and 
broken-line fit to each redshift bin individually. The break near k ~ 0.02 s km~^ is due to the change 
of regime in the splicing technique. 


'The ‘exact’ simulation with 2048^ particles per species and (100 Mpc)" box volume required 400 khrs split among 
8192 cpus for the hydrodynamical part, and a total of 300 khrs split over 128 cpus each with 4 Gb of memory for the 
extraction of the 100 thousand lines of sight. 
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2.2.3 Model of IGM temperature 

In Paper I, the relation between the density and the temperature of the IGM was derived directly 
from the simulations described in [7]. Two variables Tq and y parametrize this relation according to 
T = Tq ■ where 6 is the normalized density pi (p). For a given simulation, we measured Tq and 
y by building, for each redshift z, the T - p diagram. In the low density and low temperature region 
where the IGM lies, a linear fit of ln(r) as a function of ln(d) allowed us to determine Tq{z) and y(z). 
The parameter y(z) is monotonically and smoothly decreasing with redshift, whereas To(z) presents 
two regimes with a break at z ~ 3. The latter distribution follows notably well the measurements 
of [21]. As a consequence, in Paper I, we fixed fhe redshiff dependences fo fhaf measured in fhe 
simulations, and we lef free fwo global paramefers Tq and y. 

In fhis new analysis, we release fhe shapes of Tq{z) and y(z). Their evolufion wifh redshiff is 
modeled wifh power laws. For y we infroduce fwo paramefers: y(z - 3), fhe value of y af z - 3, 
and fhe exponenf rf of fhe relation y(z) = y(z - 3) • [(1 + z)/4]''’^. For fhe femperafure, we use 
fhree paramefers: Tq (z = 3), fhe value of fhe femperafure af z - 3, and fwo slopes (z < 3) and 
7/^“ (z > 3) fo fake info accounf fhe break af redshiff z ~ 3. 

In summary five paramefers, Tq (z = 3), (z < 3), jf° (z > 3), y (z = 3) and rf are now 
floafed in fhe fif fo model fhe densify and fhe femperafure of fhe IGM. 

2.2.4 Reionization history 

The hydrogen reionizafion history will alter the pressure smoothing scale of gas in the IGM, particu¬ 
larly at redshifts approaching the tail-end of the reionizafion at Zre ~ 6 [22]. Increasing the redshift of 
reionization allows more time for pressure to suppress small-scale structures. In the range of scales 
probed by the BOSS Lytr data. Fig. 13 of [23] shows that an increase of the redshift of reionization 
from Zre = 7 to 17 suppresses the Lya flux power spectrum at the largest modes {k ~ 0.02 s.km“^) 
by about 1% at z - 2.1 and 5% at z = 4.0. This range for Zre is extreme considering recents limits 
from Planck temperature plus low-^ polarization data, which give Zre = g (68% CL). Current 
uncertainties on Zre therefore translate into shifts of a few percent at most on the flux power spectrum, 
dominantly at high redshift and on small scales. 

In this work, we introduced two nuisance parameters that absorb the effects of different reion¬ 
ization histories. First, we have a nuisance parameter that accounts for the error caused by the limited 
simulation resolution: the slope &ik > kp of the splicing correction described in Sec. 2.2.2. As was 
already noted by McDonald et al. [23], such a nuisance parameter plays a dual role, since it also 
absorbs uncertainty related to the reionization optical depth (or, equivalently, to the smoothing of 
the baryon field). The allowed variation range of the slope (±8% at Icr for the largest modes) fully 
encompasses the maximum range of variation caused by different Zre- Secondly, we also introduced 
additional freedom compared to [23] by including a redshift-dependence of the nuisance parameter 
accounting for uncertainty in the spectrograph resolution (cf. Sec. 2.2.5). This correction is also only 
impacting small scales, with a range of about 1% at 1 cr for the redshift coverage of the BOSS data. 

An additional nuisance parameter accounting for different reionization histories would be de¬ 
generate with the nuisance parameters dedicated to the slope of the splicing correction and to the 
redshift-dependence of the spectrograph resolution. 

2.2.5 Technical nuisance parameters 

In Paper I, we identified several sources of systematics, and we tested their impact on the cosmo¬ 
logical and neutrino mass bounds we derived. In the present work, these systematics are all in¬ 
cluded as nuisance parameters free to vary in the fit. On the technical side, we varied a parame¬ 
ter to account for a possible redshift-dependence of the correction to the spectrograph resolution. 


such that the resulting correction to the spectrograph resolution is now the multiplicative factor 
Creso(^) = exp[-^^ • (ttreso + Aeso {z - 3))]. Both ffreso and ySreso are allowed to vary around 0 with a 
Gaussian constraint cr = (5 km s“^)^. 

2.2.6 Simulation nuisance parameters 

We include parameters to allow for additional freedom in our model of the IGM. We follow the study 
described in section 5.2 {impact of simplifying hypotheses) of Paper I, which we refer to for details 
on the analytical form of each correction. The multiplicative corrections are normalized to 1 for 
an impact identical to the systematic effect of Paper I. Since the cases we had assumed were quite 
extreme, we here allow the corresponding parameter to vary around 0 with a Gaussian constraint 
cr - 0.3. 

Although we remove all quasars flagged as having damped Lyman alpha (DLA) or detectable 
Lyman limit systems in their forest, we introduce a multiplicative factor 1 - [1/(15000.0 k - 8.9) + 
0.018] • 0.2 • ffDLA to account for a possible remaining contribution of high-density absorbers in the 
quasar spectra. This form is motivated by the study led by McDonald et al. [23] , and aoLA is free to 
vary in the fit. 

We model AGN and SN feedbacks by multiplicative factors of the form Cfeedback,!(^) - i^fz) + 
fSiiz) k) ■ afeedback,/, where i stands for either AGN or SN. The coefficients afz) and /3,(z) are derived 
from [24]. Over our range of interest, 0.001 < k < 0.02, the slope of the power spectrum can be 
reduced by 6% at most at z < 2.5 and 2% for z > 3.5 when adding AGN feedback, while it is 
increased by 5% and 1% for the same redshift bins respectively in the case of SN feedback. 

Fluctuations in the intensity of the ionizing background (also referred to as UV fluctuations) 
are accounted for by an additive correction proportional to the measured power spectrum at the pivot 
wavenumber kp = 0.009 s km“^, as motivated by the study of Gontcho et al. [25]. The normalization 
factor is a floating parameter. The correction is k-independent but evolves with redshift proportionally 
to the power spectrum. 

2.3 Impact of neutrino mass hierarchy 

Cosmology is often said to be mainly sensitive to neutrino masses through their contribution to the 
energy density, i.e., through Qy, and therefore only through their total mass X While this is true 
to very good approximation for the CMB spectrum, the impact of individual neutrino masses on large 
scale structures is a little more subtle. Although the main effect remains that of the total mass, matter 
power spectrum measurements have some sensitivity to individual masses because of two effects: (i) 
the detailed evolution of the background density close to the time of the non-relativistic transition of 
each species depends on individual masses, and (ii) so does the free-streaming scale of each species. 
When the individual masses are varied for the same total mass, these two effects lead respectively to 
different amplitude and shape in the small-scale matter power spectrum [4, 26]. 

Throughout this paper, we assume the three neutrino species to share a common mass equal 
to Yj wIv/ 3. Neutrino oscillation measurements, however, have shown that the three neutrino species 
have slightly different masses. According to the compilation of [27] from the combination of atmo¬ 
spheric, solar, reactor and accelerator neutrino experiments, the masses verify 6m^ - 7.54 ± 0.24 x 
10“^ eV^ and Ani^ = 2.43 ± 0.06 x 10~^ eV^ (using the authors’ formalism). The masses can follow 
a normal hierarchy (NH) with two light states and a heavier one, in which case the minimum total 
mass is X wJv = 0.06 eV. In case of inverted hierarchy (IH), the two heavy states are split by 6m^ and 
the lighter one is separated from the other two by Am^; the minimum total mass is then 0.10 eV. 

As the bounds on Y get closer to the 0.10 eV upper limit where we can distinguish between 
normal and inverted hierarchies, it becomes increasingly critical to test the impact of mass hierarchy 
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on the derived ID flux power spectrum. We therefore ran 3 simulations with identical cosmological 
and astrophysical parameters, as well as same total neutrino mass = 0.10 eV. The difference 
between the three simulations resides in the assumed neutrino mass hierarchy: the first one is pro¬ 
duced with degenerate neutrino masses (DM), the second one with normal hierarchy (NH), and the 
third one with inverted hierarchy (IH). In the latter two cases, the mass of each species is determined 
according to the squared mass differences of [27]. The masses of the individual species in each case 
are given in Tab. 3. 

Table 3: Masses of each neutrino species (in eV) for the three configurations of mass hierarchy 
considered in the text. 


Hierarchy 

mi 

m2 


Ymy 

Degenerate 

0.033 

0.033 

0.033 

0.100 

Normal 

0.022 

0.024 

0.055 

0.100 

Inverted 

0.0007 

0.049 

0.050 

0.100 


As shown in Figure 5 (left plot), we observe differences in the linear 3D matter power spectrum 
at the level of 0.1% between normal and degenerate hierarchies, and of order 0.3% between inverted 
and degenerate hierarchies, quasi independently of the scale k. The larger difference in the case 
of inverted hierarchy is explained by the fact that the total mass is essentially shared among two 
neutrinos instead of three for either the DM or the NH scenario. In the IH case, the lightest neutrino 
remains relativistic until late times (it becomes non-relativistic around znr ~ 3 for X - 0-1 eV). 
It thus contributes longer to the background density, thereby leading to a slower growth of cold 
dark matter perturbations and to a stronger overall suppression of power. In contrast, the NH scenario 
having three neutrinos of almost equal masses is closer to the case of three degenerate-mass neutrinos. 
The slope that appears at small k is an excess of power for IH compared to NH or DM (and excess 
of power of NH compared to DM), due to the presence of two (or one, respectively) higher-mass 
neutrinos, causing an earlier non-relativistic transition and thus a free-streaming damping restricted 
to smaller scales. Relatively to the degenerate case, normal hierarchy (and even more so inverted 
hierarchy) therefore exhibits an excess of power near ^nr ~ 10“^ /j.Mpc“'. The tail of this peak 
is the cause of the slope near k = 10“^ /j.Mpc“^ in Fig. 5. The effects of individual masses are 
nevertheless small, of 0.3% at most. 

The differences between hierarchies are even smaller when comparing ID flux power spectra. 
As shown in the right plot of Figure 5, the ratio of the ID flux power spectrum measured for inverted 
or normal hierarchy to the power spectrum for degenerate masses is compatible with 1 at better than 
0.05%, more than 10 times below the level of the statistical uncertainty in the simulations. Both plots 
are zoomed on the scales that Lya data can probe. The conversion between wavenumbers expressed 
in s km“' and in h Mpc“^ is redshift dependent and given by the factor H{z)/{1 + z). For the central 
cosmology used in our simulation grid and at redshifts z = 30, 4, 3, 2 and 0, the conversion factor 
equals 311, 126, 113 and 100 km s“^ Mpc“\ respectively. 

This test clearly justifies the simplifying hypothesis of degenerate masses used in this work. 
The effects of individual neutrino masses are too small to be measured with current experiments. 

2.4 Interpretation methodology 

In Paper I, we have simultaneously used Markov-Chain Monte-Carlo (MCMC) simulations with 
Bayesian inference and a frequentist interpretation. The results obtained with either method were in 
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Figure 5: Ratio of the 3D matter power spectra (left) and of the ID flux power spectra (right) mea¬ 
sured with normal (blue) or inverted (red) neutrino mass hierarchy to the power spectrum obtained 
assuming three degenerate-mass neutrinos, for a total neutrino mass f^my = Q.\ eV and identical 
astrophysical and cosmological parameters. The shaded area in the right plot illustrates the range of 
simulation statistical uncertainties for 35 k-modes over the range 0.001 - 0.020 s km“' as for data 
(cf Fig. 2). 


remarkably good agreement. In this paper, due to the increase in the number of parameters, all the 
results are derived with the faster frequentist approach (see [28, 29] ). 

Our determination of the coverage intervals of unknown cosmological parameters is based on 
the ‘classical’ confidence level method originally defined by Neyman [30]. We start with the likeli¬ 
hood £,{x, (Tx, 0), for a given cosmological model defined by the n cosmological, astrophysical and 
nuisance parameters 0 = {9\,..., 0„), and for data measurements x with Gaussian experimental er¬ 
rors (Tx. In the rest of this paper, we adopt notation, which means that the following quantity is 
minimized: 


X^{x,o-x-,&) = -2ln{£{x,o-x-,&)). (2.1) 

We first determine the minimum of ;y^(x, cr^t; 0) leaving all the cosmological parameters 
free. Then, to set a confidence level (CL) on any individual cosmological parameter we scan the 
variable 0,-: for each fixed value of 9i, we minimize again;y^(x, cTx', 0) but with n - I free parameters. 
The x^ difference, A^^(0,), between the new minimum and Xq, allows us to compute the CL on the 
variable, assuming that the experimental errors are Gaussian, 


CLi9i) - 1 - 



fx^(t',Ndof)dt, 


( 2 . 2 ) 


with 


fx^ihNdof) 


^/2^r{Ndof/2) 


(2.3) 


where T is the Gamma function and the number of degrees of freedom Ndof is equal to 1. This pro¬ 
filing method can be easily extended to two variables. In this case, the minimizations are performed 
for n - 2 free parameters and the confidence level CL(9,-, 9j) is derived from Eq. 2.2 with Ndof - 2. 
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In this paper we also combine the derived from the Lya likelihood with that of Planck. 
However, we do not directly use the Planck likelihoods. Instead, we use the central values and the 
covariance matrices available in the official Planck ^ repositories for the cosmological parameters (erg. 
Us ,0.m, Hq, Yj niv, dris/dlnk, r and r). For each parameter, we assume a Gaussian CMB likelihood 
with asymmetric Icr errors that we estimate on either side of the central value from the Icr lower 
and upper limits, thus accounting for asymmetric contours. We validated this strategy for some 
configurations by comparing it to the MCMC approach using the full likelihood. The results obtained 
for the two methods are quite similar, as already observed in Paper I. 

3 Cosmological implications of Lycr data 
3.1 ACDMv cosmology from Lya data alone 

The methodology upgrades presented in section 2 impact the constraints we derive on the parameters 
of the base ACDMv cosmology compared to Paper I, even for Lya data alone. In particular, the 
additional freedom in the model of the IGM temperature allows us to fit all 12 available redshift 
bins, thus covering 2.1 < z < 4.5, whereas we limited the analysis to the first ten bins only (2.1 < 
z < 4.1) in Paper I. As explained in the previous section, we have also more accurately determined 
the contribution of sample variance to the simulation uncertainties. Finally, we have updated the 
likelihood with more parameters to describe the IGM temperature and its evolution with redshift, 
to include a possible contribution from AGN or SN feedbacks though a parameterized correction to 
the power spectrum, and to better model the correction for the use of the splicing technique or the 
correction for uncertainty on the model of the spectrograph resolution. 

As in Paper I, we impose a Gaussian constraint on the Hubble constant Hq. This is justified by 
fhe very weak dependence of Lya dafa on Hq. Indeed, leffing Hq free in the fit does not significantly 
alter the best-fit value of any of the other parameters, but only constrains the Hubble constant to 
Hq < 77.0 km s Mpc“^ at 95% CL. In table 4, we summarize the results obtained on the fit to the 
Lya ID flux power speefrum for fhree differenl Hq consfrainfs: in column (1) we impose fhe same 
consfrainf as in Paper I (coming from Planck 2013 measuremenfs [15]), in column (2) we use a 
consfrainf faken from fhe Planck 2015 TT -i- lowP resulfs [5] and in column (3) we use a consfrainf 
faken from fhe reanalysis of Cepheid dafa by Efslalhiou [31]. There is no nofable difference in 
the values of any of the parameters whatever the constraint imposed on Hq. The resulting is 
not affected either (A;y^ <s; 1) by this constraint. The 2D constraints in the ns - (t%, Yj^v - 
and Y - fs planes are shown as the red contours in Fig. 6. The neutrino mass is correlated to 
CTg (-48%), Us (48%) and (52%). Correlations between all other cosmological parameters have 
smaller amplitudes. Astrophysical parameters, in contrast, exhibit higher correlations, with ~ 70% 
correlation between Tq, y and rf. 

The fitted values of the nuisance parameters are all well within the expected range. The correc¬ 
tion to the power spectrum due to the use of the splicing technique is fitted with a small negative slope 
utk > kp and percent-level offset at mode kp, in good agreement with Fig. 4. Both values are com¬ 
patible with 0 at 2cr.The best-fit value of the correction to the estimate of the spectrograph resolution 
indicates an overestimate ranging from 3.5 to 4.0 km s“' from lowest (mean z = 2.2) to highest (mean 
z = 4.4) redshift bin, compatible within error bars with no correction. The IGM nuisance parameters 
are all also compatible with 0, with best-fit values indicating 4%, 7% and 15% of the maximum DLA, 
SN feedback and UV fluctuation effects described in Sec. 2.2.6. There is no indication of any con¬ 
tribution from AGN feedback. The IGM temperature parameters have large error bars and are thus 

^http://wiki.cosmos.esa.int/planckpla2015/index.php/Main_Page 
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Table 4: Best-fit value and 68% confidence levels (95% upper bound in the case ofY^niy) of the 
parameters of the model fitted to the flux power spectrum P{ki,Zj) measured with the BOSS Lya data. 
The Gaussian constraint on Hq is defined by Planck 2013 cosmological results (column 1), Planck 
2015 results (column 2) or Efstathiou’s reanalysis ofCepheid data (column 3). 


Parameter 

(1) Lya + //Uaussian 

(2) Lya + 

(3) Lya + //Gaussian 


(Ho = 67.4± 1.4) 

(Ho = 67.3 ± I.O) 

(Ho = 70.6 ± 3.3) 

ers 

0.830 + 0.032 

0.831 +0.031 

0.831 +0.032 

ns 

0.939 ± 0.010 

0.939 + 0.010 

0.932 + 0.010 


0.293 + 0.013 

0.293 + 0.014 

0.292 + 0.014 

Hq (km Mpe^*) 

67.4+ 1.2 

67.3+1.0 

70.2 + 1.0 

Z niy (eV) 

<1.1 (95% CL) 

<1.1 (95% CL) 

<1.1 (95% CL) 

/Silll 

0.0059 + 0.0004 

0.0059 + 0.0004 

0.0058 + 0.0004 

/siii 

0.0007 + 0.0005 

0.0007 + 0.0005 

0.0007 + 0.0004 

To (z ^ 3) (K) 

9500 + 3500 

8900 

8500 

(z < 3) 

-2.9 + 0.4 

-2.9 + 0.5 

-3.1 +0.2 

r]^° (z > 3) 

-4.4+ 1.0 

-4.4+ 1.1 

-4.5 + 0.5 

7 

1.0 + 0.2 

0.9 + 0.1 

1.0 + 0.2 

rf 

1.0 + 0.2 

0.8 + 0.5 

0.8 + 0.2 


0.0026 + 0.0001 

0.0025 + 0.0001 

0.0025 + 0.0001 


3.734 + 0.015 

3.73 + 0.02 

3.729 + 0.005 

reduced 

0.99 

0.99 

0.99 


poorly constrained by this data set. Their values are within 1 - 2 cr of typical measurements (see, e.g. 
[21]). As explained in Paper I, there is also one additive noise correction parameter per redshift bin. 
Its value ranges from -9% to +16% with median at -2.9%, and the extreme values are obtained for 
the redshift bins at mean redshifts of 4.0 and 3.8, where the noise level is only ~ 10% at most (e.g., 
for the smallest scales) of the Lya power spectrum. 

The values of the cosmological parameters shown in table 4 are consistent within Icr with their 
value in Paper I. This indicates that the changes in methodology and data set have not modified 
our main conclusions on ACDMv cosmology. We have included additional ’nuisance’ parameters, 
which in principle should lead to increased uncertainties on fit parameters, but in most cases these 
are compensated for by the additional lever arm in redshift due to the extra redshift bins z = 4.2 
and z = 4.4 that are now included. As a result, the uncertainties on the parameters are similar or 
even slightly reduced compared to Paper I. The upper bound on the sum of the neutrino masses is 
unchanged at 1.1 eV (95% CL). The best-fit value of the sum of the neutrino masses has moved 
up slightly to 0.41 eV, although still compatible with 0 at about Icr. The main changes occur on 
erg (0.3cr downward shift) and n^ (l.Ocr upward shift). They are mostly due to our improved model 
of the splicing correction, where the broken slope dependence affects the amplitude (erg) and scale 
dependence (nfi of the fluctuations. 
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3.2 ACDMv cosmology from Lya data and other probes 

In this section, we combine the Lya likelihood (imposing no constraint on Hq) with the likelihood 
of Planck 2015 data that we derive from the central values and covariance matrices available in the 
official 2015 Planck repository. As in the previous section, we focus on the base ACDMv model, 
and we derive constraints on erg, n„„ Hq and X Results are shown in table. 5. Column (1) 
is the same as column (2) of table 4 and recalls the results for Lya alone. Column (2) is for the 
combined set of Lya and the base configuration we chose for Planck data, i.e. TT+lowP (cf. details 
in Sec. 2.1.2). The last two columns (columns 3^) include BAO data in addition, and in column 
(4) we extend the CMB measurements to TT+TE+EE+lowP. We illustrate the main 2D contours on 
cosmological parameters in figures 6 and 7. 

Table 5: Best-fit value and 68% confidence levels of the cosmological parameters of the model fitted 
to the flux power spectrum P{ki,Zj) measured with the BOSS Lya data combined with several other 
data sets. 


Parameter 

(l)Eya 

_i_ ijGaussian 

+ H^ 

{Ho = 61.3 ± 1.0) 

(2) Eya 

± Planck TT+iowP 

(3) Eya 

± Planck TT+iowP 

±BAO 

(4) Eya 

± Planck TT+TE+EE+lowP 

±BAO 

0-8 

0.831 ±0.031 

0.833 ±0.011 

0.845 ± 0.010 

0.842 ± 0.014 

ns 

0.938 ± 0.010 

0.960 ± 0.005 

0.959 ± 0.004 

0.960 ± 0.004 


0.293 ± 0.014 

0.302 ± 0.014 

0.311 ±0.014 

0.311 ±0.007 

Hq (kms^' Mpe^*) 

67.3 ± 1.0 

68.1 ±0.9 

67.7 ± 1.1 

67.7 ±0.6 

Y^niy (eV) 

<1.1 (95% CL) 

<0.12 (95% CL) 

<0.13 (95% CL) 

<0.12 (95% CL) 

Reduced 

0.99 

1.04 

1.05 

1.05 



Figure 6: 2D confidence level contours for the {(T^,ns) , {Q.m,Yj^v) (^nd (erg, ^ my) cosmological 
parameters. The 68% and 95% confidence contours are obtained with different combinations of the 
BOSS Lya data presented in section 3.1 of the Gaussian constraint Hq = 61.4 ±1.4 km s“^ Mpc“' 
and of Planck 2015 data (TT+lowP). 


The main point to note is the excellent agreement between the results derived from the combina¬ 
tion of Eya data with different sets of CMB and BAO data (columns 2-A). The consistency between 


- 14 - 
































Figure 7: 2D confidence level contours for the {cr%,ns) , and (erg, ^ my) cosmological 

parameters. The 68% and 95% confidence contours are obtained for four combinations, starting 
with Planck 2015 data (TT+lowP) alone, then adding the BOSS Lya data, the BAO data and finally 
the high-£ polarization spectra for Planck (TE and EE). 


these additional data was already discussed in the Planck cosmology paper [5]. The agreement with 
the values obtained from Lya data alone is also remarkable, and significantly improved compared to 
the previous Lya analysis (Paper I): the downward shift of erg and upward shift of noted in the 
previous section bring these two parameters closer to the value measured from CMB data. The upper 
bound on total neutrino mass is tighter than in Paper 1. We now constrain X my to be less than 0.12 eV 
(at 95% C.L.) from Lya and Planck TT+lowP, instead of 0.15 eV in Paper 1. The improvement of 
our limit is largely dominated by our improved calibration of the hydrodynamical simulations, and in 
particular by our new model of the bias induced by the splicing technique. The use of two additional 
redshift bins, and the update to Planck 2015 data also contribute to this improvement, although to a 
lesser extent. The main restriction to further improvement comes from the additional flexibility we al¬ 
lowed in the parameters describing the temperature evolution of the IGM. Our new constraint is much 
closer to the inverted-hierarchy lower bound of 0.10 eV than current CMB-based limits. For compar¬ 
ison, Planck (TT-i-lowP) alone constrains the sum of the neutrino masses < 0.72 eV, Planck 

(TT-i-lowP) -I- BAO to Yjtny < 0.21 eV, and Planck (TT-i-TE-i-EE-i-lowP) -i- BAO to ^my < 0.17 eV, 
all at 95% CE. 

The small tension on that was discussed in Paper I is still present, at the 2.3 cr level, although 
no longer affected by the presence or not of additional BAO data. We will study in Sec. 3.4 possible 
extensions to the base ACDMv models to incorporate this tension, and evaluate its impact on inflation 
models and a possible running of the scalar spectral index. The tension on has little effect on the 
constraint on X wiy because of the mild correlation between these two parameters (48% in Eya, -45% 
in Planck TT-i-lowP). As in Paper I, X is mostly correlated to erg (-48% in Eya, -95% in Planck 
TT-i-lowP), and to n„, (52% in Eya, 92% in Planck TT-i-lowP). 

It should be noted that the combination of CMB and Eya is a very efficient way of constraining 
cosmological parameters, especially ^niy. As one can see in Eigure 5, the Planck and hya+Ho 
contours in the X my - O,,, and ^ my - erg planes are complementary. The Eya data constrain Q„j 
and erg largely independently of ^ niy because they have different impact on the shape of the power 
spectrum (see discussion in §5.1 of Paper I). Eor the Planck constraints, high corresponds 
to low erg because of the suppression of power on small scales by neutrino free streaming. The 
positive correlation between and ^ my is more subtle: with Q.ch^ and Q.[,h^ well constrained by 
the acoustic peaks, raising ^ my increases the matter density at low redshift after neutrinos become 


- 15 - 























non-relativistic, and within ACDM this requires a decrease in h to maintain the well determined 
angular diameter distance to last scattering, and this in turn corresponds to higher (see, e.g., §6.4 
of [5]). The end result is that the Lya and Planck contours intersect only near ^ nty = 0. As shown 
in Figure 7, adding polarization or BAO to the Planck + Lya contours does not lead to significant 
further improvement of the constraints in these planes, at least within the ACDM framework. 


3.3 Reionization optical depth from Lya data and other probes 


In this section, we use the same Planck data sets as in Sec. 3.2, and we now consider constraints 
that can be set on the optical depth t. More precisely, we study the improvement, over what CMB 
alone can do, provided by the addition of Lya and BAO data sets, through correlations and tightened 
constraints on the cosmological parameters each is sensitive to. 

The optical depth is mainly constrained by the CMB through two effects [5, 15]. First, small- 
scale fluctuations in the CMB are damped by Thomson scattering from free electrons produced at 
reionization. This scattering suppresses the amplitude of the CMB acoustic peaks by Therefore 
Planck measures in the TT power spectrum the damped amplitude e~^'^As on small scales. This effect 
introduces a strong correlation between r and cr%. Secondly, r can be measured from CMB polariza¬ 
tion in the multipole range ^ = 2 - 30. So far the analysis presented by the Planck collaboration [5] 
used the low-resolution LFI 70GHz maps for small i. The t measurement is improved by the addi¬ 
tion of the polarization information, but this improvement does not break the degeneracy between r 
and cTg, as shown on figure 8. We can expecf further improvemenf once fhe measuremenfs of low 
mulfipoles wifh HFI polarizafion are released. 

We give below fhe consfrainfs on fhe optical depfh r and on fhe redshiff af reionizafion Zie in fhe 
framework of fhe concordance model ACDM for four combinations of dafa sefs. The Planck (TT-i-lowP) 
and fhe Planck (TT-tlowP) -i- BAO cases are faken from [5]; we add fhe Lya likelihood fo eifher case 
fo produce fhe lasf fwo configurafions. 


T-0.078 ±0.019 

^re — 

Q 0+1-8 
^■^-1.6 

r - 0.080 ± 0.016 

^re — 

10.21}j 

T - 0.080 ± 0.017 

Zre — 

10.i:j;6 

r- 0.083 ±0.015 

Zre — 

10.41};5 


Planck (TT -i- lowP) 

Planck (TT -i- lowP) -i- Lya 
Planck (TT -i- lowP) -i- BAO 
Planck (TT + lowP) -i- BAO -i- Lya 


The uncerfainfies on r and Zie are fighfened more significanfly by adding Lya fhan by including BAO. 
The values of Zre are slighfly higher fhan for Planck alone when combining eifher wifh BAO or/and 
Lya. The gain in precision is clearly due fo fhe correlafion befween fhe paramefers r and cr^ and fo a 
lesser exfenf fo fhe correlafion befween r and (see figure 8). 

3.4 Primordial fluctuations from Lya data and other probes 

In this section, we use the same data sets as in Sec. 3.2 but we extend the model to allow for additional 
parameters. In particular, the small tension between the values of preferred by Lya or Planck data 
motivates a combined fit allowing Ug to vary with scale. We thus introduce dn^/dlnk by using the 
corresponding Planck chains and adapting the Lya likelihood to include a running of We choose 
a pivot scale ko = 0.05 Mpc“', the same as in the analyses led by the Planck collaboration in order 
to allow direct comparisons. The scalar mode power spectrum is then parameterized by a power law 
with 

\ns—l + jdnJdlnk\n{k/ko) 
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Figure 8: 2D confidence level contours for the {cr^,T) and {Q.fn,T) . The 68% and 95% confidence 
contours are obtained for three combinations, starting with Planck 2015 data (TT+lowP) alone, then 
by adding the BAO data and finally the BOSS Lya data. 


As dnjd In k is not included in the grid, we take into account the variation with scale of n^ in the Lya 
data in the following way: we replace in the Lya likelihood by ng^iyaik) = + 1/2 • dn^/dlnk • 

ln(k/ko). With this approximation, we can fit CMB and Lya data with common and dafidlnk 
parameters. The pivot scale k^ is approximately in the middle of the logarithmic range of the scales 
probed by Planck. The Lya forests cover scales ranging from ~ 0.07 Mpc“^ to ~ 1.7 Mpc“\ with a 
pivot near k^ya ~ 0.7 Mpc“^. The constraint that we can derive on the running index dn^/dInk comes 
mostly from the different levels of the power spectra at the CMB and Lya pivot scales ko and k^ct- 

On the theoretical side, the simplest inflationary models predict that the running of the spectral 
index should be of second order in inflationary slow-roll parameters and therefore small, Idn^/d In k| ~ 
{ns - Vfi ~ 10“^ [32]. Nevertheless, it is possible to accommodate a larger scale dependence of 
by adjusting the third derivative in the inflaton potential, for instance. A negative running of order 
10“^ can also be implemented in inflation models with oscillations in the inflaton potential, using 
axion monodromy [33]. On the experimental side, recent CMB experiments have a mixed history 
of null-results and a-few-sigma detections of running of the scalar index. The final 9-year WMAP 
analysis found no evidence of running using WMAP alone, with dn^/dlnk = -0.019 ± 0.025 at 68% 
CL, while the combination of WMAP data with the first data releases from ACT and SPT found a 
negative running at nearly the 2cr level with dn^/dlnk = -0.022 ± 0.012 [34]. The ACT 3-year 
release measured dn^/dlnk = -0.003 ± 0.013 when combining with WMAP-7 [35]. A negative 
running was detected at just over 2o- by SPT, dn^/dlnk = -0.024 ± 0.011 [3]. The Planck 2015 
results, while roughly consistent with zero running of the scalar spectral index, indicate a ~ Icr 
preference for negative running, dn^/dlnk = -0.0084 ± 0.0082, with a slightly lower significance 
when adding HFI polarization or BAO (dn^/d In k = -0.0057 ± 0.0070), and on the contrary a slightly 
higher significance when tensor fluctuations are also allowed in addition to running (dn^/dlnk = 
-0.0126 + 0.0090). 

In the table below, we give the values of and of dn^/d In k that we find for the ACDM model 
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with a tensor-to-scalar ratio r = 0. The corresponding 2D contours are illustrated in figure 9, left plot. 


n, - 0.965 ± 0.007 
n, = 0.964 ± 0.005 
Us - 0.960 ± 0.004 
n, - 0.961 ±0.004 


dnjdlnk = -0.0084+ooo^2 

dn,/dlnA:--0.0178:°;0«54 

dn,/dln^ - -0.01491^50 
dn,/dln^--0.0152l0;0050 


Planck (TT ± lowP) 

Planck (TT + lowP) + Lya 
Planck (TT, TE, EE + lowP) + Eyor 
Planck (TT, TE, EE + lowP) + BAO + Eya 


Allowing a running of ns improves the fit;^^ by ~ 10 compared to the results obtained for the 
same set of data but without running. This is driven by the fact that a negative running of order 10“^ 
is favored both by Planck data alone, and by the tension on ns between Planck and Eya data sets. The 
four different combinations of Planck, Eya and BAO data shown above all indicate a preference for a 
negative running at the ~ 3cr level. When allowing for running, the best-fit value of ns is in excellent 
agreement with the value favored by Planck data alone. The best-fit value of dn^/dlnk is also in 
agreement (at about 1 cr) with the value measured from CMB data alone by the Planck collaboration, 
although the origin of this detection is to first order only dependent on the different values of in Eya 
and CMB data, and not on the value of running measured from either probe alone (see discussion in 
Sec. 3.5). The detection of running in the combined fit, however, could be the result of a coincidence 
between a ~ 1 cr effect in Planck, due to the mismatch between the high and low multipoles in the 
temperature power spectrum on the one hand, and a ~ 2.3 cr effect in Eya data, possibly coming from 
an unidentified systematic bias on the other hand. 

As the globalis clearly improved by letting dn^/d In k free, it is interesting to study the impact 
of this extra parameter on the determination of ^ niy in the base ACDMv model with running. As 
shown on the middle plot of the figure 9, the correlation between dn^/d In k and ^ niy is small. In the 
table below, we give the values of X and of dus/dlnk for the four configurations already studied. 


I.my < 0.65eV(95%CE) 
I-niy < 0.19eV(95%CE) 
Zmy <0.19eV(95%CE) 
Emy < 0.12eV(95%CE) 


dns/dlnk - - 0 . 0078+°°°*4 

dnjdlnk- -0.0178+°°°54 

dns/d\nk = - 0 . 00135 +°°°^® 
dn,/dlnk = -0.00141:°;°°^’ 


Planck (TT ± lowP) 

Planck (TT ± lowP) + Eya 
Planck (TT, TE, EE ± lowP) + Eya 
Planck (TT, TE, EE + lowP) + BAO + Eya 


We obtain an impressive improvement on the bound on ^ niy by including Eya on top of CMB data. 
However, the gain is slightly lower than when dn^/dlnk is fixed to zero. It may be an indication 
that part of the improvement obtained with Eya data is due to the small tension on the value of 
between CMB and Eya. Einally, even with dns/dlnk as an additional free parameter, we obtain the 
same constraint Yj^v < 0.12 eV (95% C.E.) when combining the three probes CMB, Eya and BAO. 

Eetting r free does not change qualitatively the constraint on dns/dlnk, as shown in the table 
below and as illustrated in the right plot of figure 9. In the combination, we also introduce the 
BICEP2/Keck Array-Planck (BKP) data set for reasons explained hereafter. 


ns = 0.967 ± 0.007 
ns = 0.966 ± 0.006 
ns = 0.966 ± 0.006 
ns = 0.967 ± 0.005 


dns/dlnk = - 0 . 0126 +°°';^^ 
dns/dlnk = - 0 . 0206+°°054 
dn,/dlnk--0.0im0;0«^5 
dns/dlnk = - 0 . 0200+°°054 


Planck (TT + lowP) 

Planck (TT + lowP) + Eya 
Planck (TT ± lowP) ± BKP 
Planck (TT + lowP) + BKP + Eya 


The CMB and Eya measurements provide a powerful probe of cosmic inflation through the two 
parameters and dn^/d In k as explained before. In addition, CMB polarization makes it possible to 
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Figure 9: Constraints on the scalar spectral index n^, the running dns/dlnk, and ^ niy. Left and 
middle: 68% and 95% confidence contours obtained for four combinations — Planck 2015 TT+lowP 
data alone, then adding BOSS Lya, high-i polarization from Planck (TE and EE) and finally BAO 
data. Left plot is for ACDM with running, middle plot for ACDMv with running. Right: the tensor- 
to-scalar ratio r is floated. The 68% and 95% confidence contours are obtained for four combina¬ 
tions — Planck 2015 TT+lowP, Planck 2015 TT-\-lowP and BOSS Lya, Planck 2015 TT+lowP and 
B1CEP2/Keck Array and finally all the above together. 


constrain the tensor-to-scalar ratio, r, which is directly related to the inflation field. Even if Lya data 
alone cannot measure r, it can be used in combination with CMB to improve the uncertainty on r 
over what CMB alone can do, thanks to the correlations of r with the other cosmological parameters 
and the tightened constraints Lya provides on the latter. The improvement appears on the left plot of 
figure 10 (blue and red curves). 

The uncertainty on r can also be reduced by direct measurement of the large-scale B-modes in 
CMB polarization. We therefore include the BICEP2/Keck-Array-Planck (BKP) joint analysis [16]. 
The gain provided by the B-modes for CMB data alone is visible on the left plot of figure 10 by 
comparing the blue and the purple curves. As shown in the right plot of Pig. 9, however, the BKP 
data has no impact on the measurement of dnfid In k. The addition of Lya yields further improvement, 
illustrated as the green curve in both plots. We obtain an upper limit r < 0.098 at 95% CL, letting the 
running dn^/d In k free. 

The parameters of the scalar and tensor power spectra may be estimated in the framework of 
slow-roll inflation (see [36-38]) from the value of the Hubble parameter and the hierarchy of its 
time derivatives. In the context of the Hubble flow-functions (HPP), these parameters are defined 
as: s\ - and £,+i = -Ejl(Hsi) with i > 1. The scalar spectral index the tensor spectral 

index nt and the running of the scalar spectral index dnjdlnk can be related to the three slow-roll 
parameters ei, £2 and £3 by the following equations: 

ns - I - -2s\ - £2 - 2£j - (2C -I- 3)£i£2 - C£ 2 £ 3 , (3.2) 

nt - -2s\ - 2£j - 2(C -I- 1)£i£2, (3.3) 

dn^/dln^ = -2 £i£2 - £2^3, (3.4) 

where C ^ -0.7296. The tensor-to-scalar power ratio r can be derived from the tensor index nt 
through the consistency relation nt = -r{2 - r /8 - ns)IS, obtained at second order when inflation is 
driven by a single slow-rolling scalar field. 

The middle plot of figure 10 shows the 2D contours for two slow-roll parameters £1 and £ 2 , let¬ 
ting free the running index £ 3 . This approach is equivalent, in frequentist interpretation, to marginal- 
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izing over £3 in Bayesian analysis. The constraints on these two slow-roll parameters are given below 
for four different data sets: 

£i < 0.0126 (95%CL) 

£2 - 0.0343_q QJQ 4 

Planck (TT + lowP) 

£i < 0.0108 (95%CL) 

£2 = 0.0393:«;00^^ 

Planck (TT + lowP) + Lya 

£i < 0.0067 (95%CL) 

£2 = 0.0365!«;0095 

Planck (TT + lowP) + BKP 

£i < 0.0063 (95%CL) 

£2 = 0.043i:«;°°g 

Planck (TT + lowP) + BKP + Lya 


The third parameter in the HFF formalism, £ 3 , is constrained by the measurement of dn^/d In k. 
As shown on the right plot of figure 10, Lycr data improves the determination of £ 3 . Its positive 
best-fit value directly reflects the non-zero detection of running. We should stress that a large running 
leads to values of £3 much larger than those of £ 1 , £ 2 , suggesting that the HFF expansion might not 
be convergent, and that the use of slow-roll expressions truncated at order three like in Eq. (3.4) is not 
necessarily appropriate in this context. In a future study, we will develop an analysis similar to the 
one done in [37, 39, 40] to directly measure the shape of the inflation potential that provides a good 
fit to the data, without requiring any kind of slow-roll expansion. 




0.016 

oT 

0.014 


Planck (TT+lowP) 

Planck (TT+lowP) + Ly-a 
Planck (TT+lowP) + BKP 
Planck (TT+lowP) + Ly-a + BKP 



Figure 10: Left plot: constraints on the two parameters describing primordial fluctuations (ns,r). 
Middle plot: constraints on the two slow-roll parameters (£ 2 ,£ij. Right plot: constraints on the 
two slow-roll parameters (£ 3 ,£i). The 2D 68% and 95% confidence contours are obtained for four 
combinations: Planck 2015 data (TT-t-lowP), Planck 2015 data (TT-t-lowP) and the BOSS Lya data , 
Planck 2015 data (TT-t-lowP) and the B1CEP2/Keck Array data (BKP) and finally all the data. 


3.5 Discussion on ^ my and dnftdXnk 

The joint Lya + CMB analysis presented in this paper allows us to set stringent constraints on cos¬ 
mology, providing significant improvements upon CMB alone from Planck Collaboration et al. [5] 
on two main fronts: the constraint on X my is much tighter, and a running of is measured at more 
than 3 cr. We here discuss these two results in terms of their robustness and correlations. 

The constraint on ^ my, on the one hand, comes from the measurement of cr$ in Lya data, and 
from the correlation between erg and ^ my provided by CMB. The value of erg is derived from the 
normalization of the ID flux power spectrum, and thus from its measurement at zeroth order. 

The value of dnjd In k, on the other hand, is derived from the different values of determined 
by CMB and Lya (e.g. from the measurement of the ID flux power spectrum at first order to access 


-20- 

































the slope information that determines n^), and from variations of fig within the scales probed by either 
probe (e.g. from the measurement of the ID flux power spectrum at second order). It is therefore 
more sensitive to systematic effects in the measurement of the flux power spectrum. For instance, 
the slope of the ID flux power spectrum is sensitive to the modeling of instrumental effects (such 
as spectrograph resolution) as well as of physical contributions affecting the intergalactic medium 
(such as SN or AGN feedback, UV fluctuations, etc.). These effects have been included in this work 
through nuisance parameters that are fitted along with the relevant cosmological parameters. This is 
nevertheless a delicate task, not free of any possible yet-unaccounted-for additional systematic which 
could affect the determination of 

We can estimate the value of dn^/d In k from the difference of scale factors at the CMB and Lya 
pivot scales, kcMB - 0.05 Mpc“^ and kiya ~ 0.7 Mpc“^ respectively. Given the definition of running 
of Eq. 3.1, and the values of determined separately from CMB or Lycr data, we estimate dn^/d In k 
to be approximately -0.02. This is in agreement with the best-fit value -O.OlVS^QQQ^g (cf. Sec. 3.4), 
thus confirming that running is indeed detected in this work mostly from the different levels of in 
CMB and Lya data, and thus mostly from a first order measurement of the Lya power spectrum. Any 
unidentified systematic uncertainty that would resolve the tension on Ug would thus simultaneously 
annihilate our detection of dn^/d In k. 

Linally, let us study the impact on X of the 2.3 cr tension on Ug. In Paper I (Sec. 4.2.2), we 
showed that it did not affect the result since the tension on Ug is roughly the same for small (0.1 eV) or 
large (0.3 eV) neutrino masses. In the present work, we infer a similar conclusion from the fact that 
Hg and Yj niy have little correlation, both in Lya and in CMB data (cf. Sec. 3.1 and 3.2). As a final test, 
we implemented a dedicated MCMC for CMB data, allowing both ^ niy and dng/dlnk to vary. The 
result, given in Sec. 3.4, indicates a null correlation between these two parameters. Errors or limits 
on both parameters are also in excellent agreement with the values obtained when only one of them 
at a time is included in the fit, confirming the absence of degeneracy between ^ niy and dng/dlnk. In 
conclusion, if in the future a systematic effect is found that affects the determination of Ug, it would 
directly alter the determination of dug/dlnk but would not modify significantly the limit on X "Jy 

4 Conclusions 

In this paper, we present an update on the constraints we derive on several cosmological parameters 
using Lya data, either taken alone or in combination with CMB and BAO data. We improve upon our 
previous study of Paper I [2] on several fronts that we summarize below. 

We improved the likelihood describing the Lya data in several ways. We relaxed our model of 
the IGM, now modeling the redshift-dependence of the IGM temperature Tq and its 6 dependence 
(i.e., the logarithmic slope y) as, respectively, a broken and a single power law. The likelihood now 
includes additional freedom to account for the systematic uncertainties that we identified. We ran 
new hydrodynamical simulations that allow us to improve our model of the impact of the splicing 
technique on the ID power spectrum. This splicing technique is used in the simulation grid to mimic 
large high-resolution simulations equivalent to 3092^ particles per species in a (100 Mpc)^ box. 
We also better estimate the contribution of sample variance to the simulation uncertainties. Using 
this updated likelihood, we set a robust upper bound on the sum of the neutrino masses Yj^v < 
1.1 eV (95% C.L.) from Lya data alone, prioritizing over all known systematics. Compared to the 
Planck 2015 measurement [5], the constraints we derive on cosmological parameters and erg are 
in excellent agreement, but we note a 2.3 cr discrepancy on Ug. 

In a second step, we combine our new Lya likelihood with the MCMC chains from Planck 
2015. In the context of a flat ACDMv cosmology, this combined study leads to the tightest limit 
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published on with an upper bound of 0.12 eV (95% C.L.) using Lya and Planck TT+lowP 
data. The addition of BAO does not further improve this limit. 

In a third step, we investigate the improvement, over CMB data alone, that Lyor data can provide 
on other parameters. In particular, we show that through the correlation between the reionization 
optical depth t and erg, Lya data can reduce the uncertainties on the measurement of r. The central 
value remains fully compatible with the one derived from CMB. 

Finally, we consider extensions to the ACDMv cosmology, in particular focusing on a possible 
running of the scalar spectral index. We note a clear improvement of the total when allowing for 
running, in agreement with the observed tension on Ug between Planck and Lya data sets. However, 
this improvement, and the subsequent 3 cr detection of running, could result from a coincidence 
between a ~ 1 cr effect in Planck (due to the mismatch between the high and low multipoles in the 
temperature power spectrum) on the one hand, and a ~ 2.3 cr effect in Lya data (possibly coming 
from an unidentified systematic effect) on the other hand. Allowing for dn^/dlnk to vary in the ht 
does not change our best limit on Yj but in that case BAO data are needed to reach that best limit. 
The combination of Lya with Planck (TT, TE, EE + lowP) and BAO still leads to ^ nty < 0.12 eV at 
95% C.E., even when allowing for running. 

As an interesting consequence of our precise measurement of dn^/dlnk, we can derive con¬ 
straints on the inflation potential. In the context of slow-roll inflation, we use the second-order 
expansion in Hubble flow-functions to relate the slow-roll parameters ei, q and 63 to the power 
spectrum parameters rit, and dnjdlnk. The uncertainties on the inflation parameters are reduced 
signihcantly compared to their measurement with CMB data. 
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